Clustering and collisions of heavy particles in random smooth flows 
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Finite-size impurities suspended in incompressible flows distribute inhomogeneously, leading to a 
drastic enhancement of collisions. A description of the dynamics in the full position-velocity phase 
space is essential to understand the underlying mechanisms, especially for polydisperse suspensions. 
These issues are here studied for particles much heavier than the fluid by means of a Lagrangian 
approach. It is shown that inertia enhances collision rates through two effects: correlation among 
particle positions induced by the carrier flow and uncorrelation between velocities due to their finite 
size. A phenomenological model yields an estimate of collision rates for particle pairs with different 
sizes. This approach is supported by numerical simulations in random flows. 
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U 



> 

o 
o 
o 



> 
X 



I. INTRODUCTION 

The dynamics of small impurities, such as droplets, 
dust or bubbles, transported by an incompressible flow is 
much more complex than that of point-hke fluid tracers. 
This is due to their finite size and to their mass density 
being different from that of the carrier fluid. As a con- 
sequence of their inertia the dynamics of such particles 
is dissipative, leading to inhomogeneities in their spa- 
tial distribution. This phenomenon, frequently referred 
to as preferential concentration, has been observed for a 
long time in experiments (see Ref. Q for a review) . Sus- 
pended particles typically interact through collisions and 
chemical or biological processes. It is therefore very im- 
portant in a large spectrum of applications to quantify 
the effects of inertia on these interactions. Le us mention 
for instance the problem of estimating the time scales of 
rain initiation in warm clouds. ^ ''■^ Other examples are 
the problem of microorganism predator-prey encounters 
in turbulent flowS j^i^i^ and the enhancement of chemi- 
cal reaction rates for active particles suspended in fluid 
flows. ^'^ Being a common characteristics of many of the 
above examples, here we focus on very dilute suspen- 
sions of very small particles that are much heavier than 
the carrier fluid. They are moved by the fluid through a 
viscous drag, whose nondimensional characteristic time, 
the Stokes number S, is proportional to the square of 
their radius. In most cases, particle sizes are below the 
smallest characteristic scale of the flow, where the fluid 
velocity is spatially smooth. We are interested in the in- 
teractions taking place at those scales and we therefore 
consider random smooth flows. 

In the last few years, much effort has been 
devoted both from a theoreticaliSiii*i2iiii4 and 
numericaliSiiSiiLi& point of view to understand and 
quantify the enhancement of collision rates induced by 
inertia. Two mechanisms have been identified: prefer- 



ential concentration increases the probability for two 
particles to be at a colliding distance)24i& detachment 
from fiuid trajectories may enhance the relative velocity 
between two approaching particles 4^*^^ Previous works 
mostly treated independently these two mechanisms. 
This is justified in the asymptotics of very small and 
very large inertia. In the former, particles are almost 
tracers so that their relative velocity is given by the 
fluid velocity shear; the discrepancy from a uniform dis- 
tribution is responsible for an enhancement of collisions. 
In the other limit particle motion is essentially ballistic: 
they distribute uniformly but they may reach the same 
position with very different velocities ~ this is known as 
the sling effect.^ 

Several attempts have been done to bridge the gap be- 
tween these two asymptotics. For instance, Kruis and 
Kustersii proposed an interpolating formula for the typ- 
ical particle relative velocity, by summing together the 
effect of shear due the fluid velocity and the acceleration 
induced by inertia. This model, though reducing to the 
result of Saffman and Turneni^i in the tracer limit and of 
Abrahamson^^ for very large inertia, does not take into 
account the effects of preferential concentration. Some 
improvements have been proposed in Refill by using two 
analytical models for the fluid-particle velocities, but still 
neglecting the effect of particle clustering. 

We propose here a new approach which treats these 
two mechanisms in a coherent manner. Of course, we 
do not to provide a complete analytical model for the 
collision rates, our aim is to introduce a consistent phe- 
nomenological framework where the quantities relevant 
to its estimation are stressed out. We mainly make use of 
tools borrowed from dynamical systemsiS to reinterpret 
inertial effects. The basic idea is to take into account the 
full position- velocity phase-space dynamics of particles. 
Preferential concentration of particles is then understood 
as the convergence of trajectories toward a dynamically 
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evolving attractor in phase space. Folding of the attrac- 
tor in the velocity direction is responsible for the increase 
of the velocity differences between particles. The statis- 
tical properties of the fractal set are determined by the 
carrier flow and depend on the particle radius through the 
Stokes number. This approach has a natural extension to 
polydisperse suspensions where particles are distributed 
over sizes: the effects of inertia on the relative motion 
of particles can there be studied in terms of correlations 
between attractors labeled by different Stokes numbers. 

As we focus on very dilute suspensions, the dynami- 
cal effects of collisions are here neglected, i.e. we work 
within a ghost particle approximationiSSiSi This allows 
us to investigate the phase space in terms of the La- 
grangian dynamics, that is following the trajectories of 
particle pairs. Thus clustering can be characterized in 
terms of the probability, Pi.2{r), for the pair separation 
to be below a distance r; in the same fashion we introduce 
the "rate of approach" , defined as the fraction of particle 
pairs at a distance r that approach each other during a 
unit time. This is given by the average of the negative 
part of the radial velocity difference at a separation r. 
As we shall see, the r-dependence of these two quantities 
sheds light on and characterizes the various regimes of 
particle motion. Moreover, the approaching rate allows 
to estimate, once r is fixed to be the sum of the particle 
radii, the collision rates. 

In the monodisperse case -Pi. 2(f) behaves as a power 
law with an exponent displaying a non-trivial dependence 
on S. This exponent is equal to the space dimension, 
d, in both the very small and very large Stokes number 
asymptotics, where particles distribute uniformly. The 
scaling behavior of the rate of approach results from the 
joint effect of clustering and velocity difference between 
particles. The latter behave very differently in the two 
asymptotics: for 5 <C 1 it is proportional to r, while for 
5 » 1 it becomes independent of the separation. For in- 
termediate values of the Stokes number, correlations be- 
tween particle separations and velocity differences leads 
to non-trivial scaling behavior for the approaching rate. 

As to polydisperse suspensions we show that for the 
relative motion of two particles with different Stokes 
numbers, a critical separation r^, oc 1 5*1 — 5*21 is singled 
out. Below it, the two motions are essentially uncorre- 
lated. Correlations, due to the fact that particles are sus- 
pended in the same flow, show up for length scales above 
r^. The origin of this characteristic length is understood 
in terms of the pair-separation dynamics, which is domi- 
nated by the Stokes difference for r < (the accelerative 
mechanism), and by the fluid velocity (the shear mech- 
anism) when r > r^,. This crossover length separates in 
two distinct regimes the scale dependence of both the 
probability distribution of particle separations and their 
rate of approach. 

Exploiting the phenomenological understanding, sup- 
plemented by numerical computation, of the r- 
dependence of the approaching rate (for equal-size and 
different-size particle pairs) we finally propose a semi- 



quantitative, phenomenological model for the effective 
collision kernel in polydisperse suspensions. 

The paper is organized as follows. In Section ^1 ^f- 
ter recalling the equations of motion, we discuss the ba- 
sic ingredients and the validity of the approach based 
on Lagrangian statistics and ghost collisions. In Sec- 
tion model flows used for numerical illustrations are 
described. In Section Hvl and IVl we investigate the scal- 
ing behavior of the pair separation probability and of 
the approaching rate in monodisperse and polydisperse 
suspensions, respectively. In Section I VII after a brief re- 
view of earlier investigations and models, we propose a 
phenomenological model for the effective collision kernel 
in polydisperse suspensions. Section IVIII is devoted to 
discussions and conclusions. 



II. DYNAMICS AND STATISTICS OF DILUTE 
SUSPENSIONS 

The dynamics of very dilute impurities suspended in an 
incompressible flow u is described by a standard model 
proposed by Maxey and Riley,^^ which was derived under 
the following assumptions. The particle radius, a, must 
be much smaller than the Kolmogorov length 77. Parti- 
cle Reynolds number has to be small enough to ensure 
that the impurity is surrounded by a Stokes flow. The 
suspension should be very dilute so that both the hy- 
drodynamical interactions between particles, and their 
feedback on the carrier fluid can be neglected. Since we 
focus on the effects induced solely by inertia, we ignore 
gravity. Here we address the problem of particles much 
heavier than the carrier fluid, so that this model reduces 
to the Newton equation)^ 



d^X IdX 1 , 



(1) 



where X{t) denotes the trajectory. The response time 
T, frequently referred to as the Stokes time, is related to 
the particle radius a by 



2 pr,a 

T = — ■ 



9 pfV ' 



(2) 



where the mass density ratio between the particle and 
the fluid Pp/pf is assumed to be very large (e.g. for wa- 
ter droplet in air Pp/pf ~ 10'^ and for aerosols Pp/pf « 
10'^— 5 X 10'^); f is the fluid kinematic viscosity. It is worth 
reminding that, in the asymptotics of very heavy parti- 
cles, the prerequisite of small particle Reynolds number 
does not restrict too much the admissible range for the 
values of the Stokes time. Indeed the nondimensional 
Stokes number S — t/t,i — 2ppa? /{9rj^ pf) (defined as the 
ratio between the particle response time and the shortest 
characteristic time scale of the turbulent fluid flow, i.e. 
the eddy turnover time associated to the Kolmogorov 
scale) may vary in a range from 0, the tracer limit, to 
about 10 for particle sizes of the order of rj/10 with large 
mass density ratio. 
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In terms of the nondimensional parameter the equa- 
tion of motion JQ) can then be rewritten as 



- = ^iuiX,t)^V), (3) 

time being now rescaled by t,,. Here, we expUcitly intro- 
duce the particle velocity V to emphasize that, at vari- 
ance with tracers, particle dynamics takes place in the 
(2 X d)-dimensional position-velocity phase space {x,v). 
It should be noted that the dynamics defined by Eq. Q 
is uniformly contracting at a rate —d/S, so that particle 
trajectories generally concentrate in phase space onto a 
dynamically evolving attractor. In the large time asymp- 
totics the phase-space density of particles, solution of the 
Liouville equation associated to ((Hjl, becomes singular 
with its support on this attractor. Its statistical prop- 
erties are usually multifractaliiSi 

We consider a collection of particles with different 
Stokes numbers embedded in a flow defined in a finite 
domain of size L. Even though the particles are car- 
ried by the same fluid flow, they converge to different 
attractors according to their Stokes number. A suitable 
characterization of suspension dynamics is given by the 
instantaneous phase-space densities fs{x,v,t). Here, fs 
is normalized by the total number of particles with Stokes 
number S, so that it can be interpreted as the probability 
density to find at time t a particle with Stokes number 
S at position x and with velocity v and for a given fluid 
flow realization. 

In experiments one often has access to the distribu- 
tion of particles in the position space only. The lat- 
ter is obtained by integrating over particle velocities, 
ns{x,t) = J dvfs{x,v,t). For tracers, the density 
no (a;, t) is uniformly distributed over the domain, so that 
the coarse-grained density no(r, t), obtained by integrat- 
ing over small volumes of size r, scales as no(r, i) ~ r**, 
being d the spatial dimension. On the contrary, for 
inertial particles, typically ns{r,t) ^ ^Di{s) -^^ji^jj jj-^ 
markedly smaller than the space dimension d. The in- 
formation dimension Di is one among the dimensions 
Dq that characterize the scaling properties of multifrac- 
tal densitiesi^^*^ For particle pairs dynamics the relevant 
quantity is the correlation dimension D2{S), which mea- 
sures the scale dependency of the probability that 
two particles on the same attractor are separated by a 
distance smaller than r. Deviations of D2{S) from the 
value d corresponding to a uniform distribution are im- 
portant to quantify the weight of clustering (induced by 
inertia) on particle-to-particle interactions. 

The radial distribution function g{r) is frequently used 
in the literature to characterize clustering of inertial 
particlesii^ This quantity can be expressed as the ra- 
tio between the actual number of particles inside an in- 
finitesimally thin shell of radius r centered on a given 
particle and the number that would be expected if the 
particles were uniformly distributed. It is easily checked 



that the radial distribution function behaves as r^^'"^)"''. 
For particles uniformly distributed, D2{S) — d, so that 
g{r) ~ 0(1). On the contrary, when D2{S) < d, the sig- 
nature of particle clustering, g(r) diverges as r — 0, as 
it was numerically demonstrated in Refs. ,16tal7i 

It is worth mentioning here that clustering is also found 
in tracers advected by compressible flows. When the 
compressibility is sufficiently large, particles can collapse 
onto point-like attractors: this is usually referred to as 
the strong compressibility regime,— In such flows, particle 
inertia may induce a further enhancement of clustering, 
leading to a strong compressibility regime even at lower 
values of the compressibility.^'' 

Suspensions generally involve local interactions among 
the particles: dry dust scatters through clastic or in- 
elastic collisions depending on their impact velocities— 
or water droplets coalesce to form rain dropsi^ In dilute 
suspensions the dominant interactions are the binary col- 
lisions taking place when two particles with radii oi and 
02 {viz. Stokes numbers Si and S2), are at a distance 
r fii ai + 02- It has been proposedi22i2i that for very di- 
lute suspensions, as considered here, the collision rates 
can be estimated by using the so-called ghost collisions 
approach. The main idea of this method is to let the 
particles overlap after they come across a distance equal 
to the sum of their radii. Neither the particle sizes, nor 
their velocities are modified after such an event: the colli- 
sions are not physically performed but just recorded and 
counted to estimate the rate at which they happen. To 
understand the relevance of such an approach let us per- 
form the following gedanken experiment: we throw inde- 
pendently a very large number of particle pairs inside a 
given domain and let them evolve in different random re- 
alizations of the carrier flow until the two particles form- 
ing these pairs collide. We define here the collision rate 
as the rate at which the number of particle pairs which 
have not collided decreases with time. The dynamics 
starts with transients during which the probability that 
two particles collide is strongly influenced, for instance, 
by the choice of the initial velocities and positions of the 
particles. When the dynamics of the particles has con- 
verged to a statistical steady state, the collision rate at- 
tains an asymptotic value which depends on the Stokes 
numbers of the particles and on the properties of the 
carrier flow. This limiting value is that measured when 
assuming ghost collisions. The latter are hence expected 
to give a good estimate of collision rates in those situa- 
tions where the time of convergence to a steady regime is 
much shorter than the mean free time between collisions 
of particles, as it is for instance the case when consid- 
ering very diluted suspensions. The framework of ghost 
collisions have already been compared to several other 
collision schemes iS2*2ii2S As one could easily expect, it 
was found to be of particular relevance when the volume 
fraction of particles is very low. Note finally that be- 
side justifying the use of this approach, the high-dilution 
hypothesis has another important aspect: it permits to 
neglect the feedback effect of particles on the fluid flow 
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(the so-called reverse coupling) and to focus on under- 
standing the enhancement of collision rates solely linked 
to particle inertia. 

In order to evaluate collision rates, we introduce the 
rate at which particles at a distance r approach each 
other, which can be written in terms of the phase-space 
densities as 

\J\xi—X2\—r I 

(4) 

where = {x^v) and the phase-space integral is per- 
formed both on velocities and on positions satisfying 
\xi — = \r\ = r, and Wj^j denotes the negative part 
of the radial velocity difference {v2 — Vi) ■ r. The aver- 
age ([■■•]) is performed over the realizations of the fluid 
velocity field. It easily understood that, for r = a\ + 02, 
the rate of approach ki^2 gives the fraction of particles of 
size ai and 02 colliding per unit time inside the domain. 
However, as we shall see, the full r-dependence of ki,2 
is very informative on the dynamics and thus interesting 
by itself. 

The role of inertia in increasing the rate Q with re- 
spect to simple tracers can be understood easily for equal- 
Stokes particles. Si = S2 — S. First, the clustering 
due to inertia enhances the probability for particles to 
be close to each other. As emphasized in Ref. 3, this 
mechanism is very important in the asymptotics 5 <C 1. 
Second, an additional increase of the rate of approach 
Ki_2 may be induced by large radial velocity differences. 
As a consequence, non-trivial dependencies of ki^2 on the 
scale separation r will emerge. This effect is known to be 
very important at large S, where caustics appear in the 
phase-space distribution of particles. Such singularities 
correspond to physical-space positions where a folding of 
the attractor leads to the presence of different branches in 
velocity at the same position, and hence, to large veloc- 
ity differences between close particles. Such a folding in 
velocity is the basic mechanism leading to what is known 
as the sling effect 

Lagrangian statistics 

Clearly, a direct (2 x d)-dimensional numerical integra- 
tion to determine phase-space densities is a difficult task. 
We thus adopt a different strategy by using a Lagrangian 
approach. Indeed, investigating two-point correlations of 
densities, as needed for binary collisions, is equivalent to 
study the relative motion of particle pairs. Moreover the 
average over velocity configurations can now be replaced 
by a time average»2^. To evaluate the rate of approach 
Ki,2, it is sufficient to follow two particles with Stokes 
numbers 5*1 and 6*2, and to compute the time average 
of the negative component of the radial velocity differ- 
ence conditioned on the pair distance, r. Numerically, 
it is more convenient to condition such averages on hav- 
ing the two particles at a distance smaller than r, rather 



than exactly equal to r. This amounts to work with cu- 
mulative quantities that are by far less noisy. We hence 
define 

/Ci,2(r) = |I^i,2| H{-Wia) H{r - |iJi.2|) , (5) 

where R12 = Xi - X2, W1.2 = {Vi - V2) ■ Hi,2/|fii,2|, 
H is the Heaviside function and the overbar denotes La- 
grangian time averaging; the rate l@J is then obtained as 
i^iA'^) = drlCi,2{r)- 

To quantify the balance between the two effects en- 
hancing collisions mentioned above, namely clustering 
and large velocity differences, we make use of the proba- 
bility density function Pi^2('') for the pair separation. It 
is clear that this quantity is related to the radial distri- 
bution function by 51, 2('') = L'^Pi,2{r)/'r'^~^ ■ Similarly, 
instead of measuring directly pi,2(?'), we consider the cu- 
mulative probability 



Pi,2(r)=i/(r-|fii,2|), (6) 

that is the probability that two particles are at a distance 
closer than r. It is clear that ^1.2(7") = drPi,2(r) (note 
that while Pi^2 is adimensional, pi^2 has the dimension 
of the inverse of a length.) 

For separations r much smaller than the smallest char- 
acteristic length scale of the fluid flow (e.g. the Kol- 
mogorov scale 77 for turbulent flows), both pi 2 and Ki^2 
are expected to display power law behaviors: 

/ r \ f^(.Si,S2) 
PiMr) - [2) , (7) 

/ r \ ■y{Si,S2) 
«i,2(r) ^ (-) . (8) 

Therefore, the corresponding cumulative quantities Pi. 2 
and /Ci_2 do also behave algebraically with exponents 
/i + 1 and 7 -I- 1. Note that these quantities are sym- 
metric under particle exchange (1,2) 1-^ (2,1), so that 
^i{Si, S2) = ^l{S2,Sl) and 7(5"!, S'2) = 7(6*2, Si). 

For equal-Stokes particles (5*1 = 5*2 = 5'), n{S) + l is ex- 
actly the correlation dimension D2{S) of the projection of 
the attractor onto the physical space. For generic pairs, 
the exponent fi{Si,S2) + 1 measures the correlation be- 
tween the distributions of particles of the two types. For 
smooth flows, the exponent 7 is constrained between fi 
and fi+l. The lower bound is attained when particle ve- 
locities and their separation r are uncorrelated, while the 
upper bound correspond to a velocity difference between 
particles proportional r. 



III. MODEL FLOWS 

Since collisions occur at distances comparable with 
the particle radii, therefore typically much smaller than 
the turbulent viscous scale 77, we can limit ourselves to 
smooth incompressible fluid velocity fields. In order to 
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perform long time averages at a reasonable cost, we con- 
sider synthetic random flows. The physical mechanisms 
are not expected to be very sensitive to the space di- 
mensionality. The only important effect is that the prob- 
ability of finding two close particles decreases with the 
dimension d. We therefore choose two-dimensional flows 
in a finite domain of size L. The validity of the results 
presented in this work should extend, at least qualita- 
tively, to more realistic velocity fields. 

First, we consider an isotropic, homogeneous and 
Gaussian flow. The Fourier modes of the fluid velocity 
are Ornstein-Uhlenbeck processes satisfying 

—7r = Mfc+CfcCfc, (9) 



where the ^^'s are independent white noises and the 
correlation time r/ is independent of the wavevector k. 
The process © is defined over 8 modes and the con- 
stants Cfc's have been suitably chosen to ensure statis- 
tical isotropy at small scales. This velocity field, being 
smooth in space and continuous in time, is expected to 
mimic the dissipative-range dynamics, when neglecting 
intermittency effects; L should then be understood as 
the Kolmogorov scale ry, and the correlation time t/ as 
the Kolmogorov time Note that the existence 

of random dynamical attractors for heavy particles, and 
hence of a statistically steady state, has been recently 
proved for this model flown^i The advantage of using a 
reduced number of modes is that the Fourier summation 
can be directly performed at particle positions, enabling 
to resolve the finest scales of the dynamics without inter- 
polation. Particles evolve according to Eq. in a peri- 
odic domain of size Lx L (here L = 2tt). Time marching 
is performed by a fourth-order Runge-Kutta scheme for 
times up to 10^ — 10^ r/. The Stokes number is defined 
a.s S = t/tj. 

Second, we investigate the alternating shear flow: 
u, = 0-u2= ^"^1(^ -^1) for t e [(n-i)T, („-i)T] , 

where periodic (of period L) boundary conditions are im- 
posed in {xi,X2), and [/„ — is chosen randomly at 
each time interval T with equal probability. In this case 
the Stokes number is defined as S — t/T. The reason 
for studying such a model flow is twofold. First, it is 
interesting for testing the robustness of the results, at 
least at a qualitative level. Second, this model flow has 
the nice feature that the integration of particle motion 
can be done explicitly in each time interval of duration 
T/2, allowing for extremely long time averages (up to 
108 _ 109 T). 



IV. LOCAL DYNAMICS OF MONODISPERSE 
SUSPENSIONS 

We start considering the asymptotics of small Stokes 
numbers, that is when the particle response time is much 
smaller than any typical time scale of the flow. A simple 
perturbative expansion of implies that inertial par- 
ticles behave as fluid tracers evolving in a slightly com- 
pressible flow whose divergence is proportional to 5*.— 
At small separations r <C i, being the carrier flow spa- 
tially smooth at these scales, the radial velocity difference 
of two particles is proportional to r. Hence, for equal- 
size particles 7 = /i -I- 1. For a vanishing Stokes number 
S = 0, the particles behave as tracers in an incompress- 
ible flow, therefore they distribute uniformly in space, 
meaning that ^ = d — 1 and "f = d. For small yet flnite 
S, if the velocity field is isotropic and homogeneous, one 
expects the discrepancy from the uniform distribution to 
behave as 5^, and in particular ^{S) ~ c? — 1 — aS^ {a 
is a fiow-dependent constant) i22i2^ In terms of the radial 
distribution function this means that g{r) diverges in the 
limit r — > as r~°'^ . This behavior has been indeed ob- 
served in Gaussian, random incompressible flows^i and in 
direct numerical simulations of three-dimensional turbu- 
lent flows. ■''^ Therefore, in the case of small Stokes num- 
bers, the main effects of inertia on particle interactions 
stems from preferential concentration, which increases 
the probability of flnding close particles Ji2£ 

The other asymptotics is also easy to understand: for 
5^1, particles have a nearly ballistic motion. They are 
thus expected to distribute uniformly in the whole do- 
main, i.e. fi ^ d—1. Moreover, even close particles have 
typical velocity differences that do not depend on their 
separation, i.e. j ^ 11 ^ d — 1. Therefore, even though in 
both large and small Stokes numbers asymptotics particle 
densities are uniform, the rates of approach are markedly 
different. 

The exponents fi and 7 computed numerically for the 
two flows described in Section lTTll are shown in Fig.QJas a 
function of the Stokes number S. At both small and large 
values of S the expected behaviors are rather clear: when 
S 1, there is coincidence between 7 — 1 and fi (velocity 
differences are proportional to the separation), while for 
S :S> I, one observes that 7^2 (velocities are almost un- 
correlated) . For intermediate Stokes numbers non-trivial 
dependencies between particle separations and velocity 
differences appear. In particular, notice that has a min- 
imum (meaning maximum of clustering) for a finite value 
of S that, as far as we know, cannot be predicted by any 
present theory. Moreover, close to this minimum 7 — 1 
noticeably deviates from fi. As we shall see in Sec. I VII 
this implies a further increase of the collision rate in this 
range of Stokes. It is also worth stressing that both mod- 
els lead to qualitatively similar results, bringing evidence 
for the robustness of these features. Note that this prop- 
erty Note that this property has already been observed 
numerically and experimentally. ^ 

The non-trivial relations between 7 and fi at varying 
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S can be understood in terms of two competing effects 
in phase space: the folding of the attractor in the in- 
direction and its tendency to contract toward the surface 
defined by the instantaneous fluid velocity. The typical 
rate of this relaxation is given by the Stokes time. For 
5 = it becomes infinitely fast, preventing folding so 
that the particle velocity field is mostly mono- valued (see 
Fig. 12^) . As the Stokes number increases the probability 
of finding particles at the same position with different 
velocities becomes larger (see Fig.l^j). Such points, once 
projected onto position space, lead to what can be seen 
as self-intersections of the attractor projection. This phe- 
nomenon was pointed out as the "sling effect" in Ref. 0. 
We stress that the numerical observations reported in 
Fig.fflshow that this phenomenon starts to be important 
already for relatively small Stokes numbers, as signaled 
by the deviation of 7 — 1 from /i. 

To summarize, as we expected, the two quantities of 
interest to measure clustering and collisions have, for par- 
ticles with the same Stokes number S, a power-law be- 



havior at small separations r: 



Ps{r) 



L \L, 



(11) 
(12) 



The constant Cs depends on S and on the statistics of 
the velocity gradients of the carrier fluid. For three- 
dimensional turbulent flows it may also depend on the 
Reynolds number. V5 is a typical velocity of the parti- 
cles with Stokes number S. For small S it is of the order 
of the root-mean-square velocity of the carrier flow, u„ns- 
For S* > 1 it can be shownS^ that Vg ~ UrmsS"^/'^. 



V. AN EXTENSION TO POLYDISPERSE 
SUSPENSIONS 

We now investigate suspensions with particles having 
different Stokes numbers. To this aim we consider the 
equations governing the relative motion of two particles 
with Stokes numbers S'l and 52. The separation R and 
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FIG. 1: (a) Scaling exponents \i of the density correlations 
and 7 of the rates of approach as a function of the Stokes 
number for the time-correlated Gaussian random flow. Time 
averages are performed over 5 x 10® r/ (b) Same as in (a) for 
the shear flow with T = L = 1 , U= -\/23. Time averages are 
over 2>d0*' T. In both cases the exponents are measured as the 
mean logarithmic derivative of K, over two to three decades in 
r. 




FIG. 2: Modulus of the particle velocities as a function of 
their positions for two different values of the Stokes number 
(a) S — 10""^ and (b) 5* = 1. At small S the surface identified 
by the particles velocities is very close to that of the modulus 
of the Eulerian velocity field, meaning that at fixed spatial 
position the distribution of the particle velocity is sharply 
peaked on the fluid one. For larger Stokes numbers, the at- 
tractor folds in the velocity direction allowing particles to be 
very close with very different velocities. 
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FIG. 3: Snapshot of the positions of A'' = 4 x 10^ particles 
associated to two different values of the Stokes number, S — 
0.4 (black) and S — 0.5 (gray, red online) for the random 
flow @. The upper- left inset shows a zoom illustrating the 
effects induced by the difference in Stokes number (see text 
for details). 



the relative velocity W evolve according to 
dR 
~dt 

dW I l^u-W 9 u-V 



dt 



T 1-6*2/4 T 1-612/4 



(13) 



where u = (tti + U2)/2, Au — Ui — U2 and V ~ (Vi + 
V2)/2. Here and in the sequel, t = (n + T2)/2, S = 
{Si + S2)/2 and 6 = - t2)/t = AS/S. 

The two terms in the right-hand side of 113|l are as- 
sociated to different effects. The first corresponds to the 
relaxation of the relative velocity to the fluid velocity dif- 
ference. The second is proportional to the difference in 
Stokes numbers, and therefore it vanishes for equal-size 
particles. The characteristic length-scale 



r. = L\e\ 



(14) 



distinguishes two ranges of scales. When R = \R\ > r*, 
the particle velocity difference is driven by the fluid ve- 
locity difference | Au| R while when i? < r^, it is driven 
by the term originating from the difference in particle 
response time. As a consequence for i? < r^, particle 
motion is uncorrelated, while above correlations in- 
duced by the fact that particles are transported by the 
same velocity field show up. This mechanism is exempli- 
fied in Fig. 121 where the simultaneous snapshots of two 
populations of particles characterized by different Stokes 
numbers are shown: at large scales the two distributions 
look essentially the same, while a closer inspection reveals 
their differences below the crossover length. 

Two cases need to be distinguished: small Stokes num- 
ber differences, corresponding to ^ L for which an in- 
termediate asymptotic range is present, and finite Stokes 



number differences for which 
termcdiatc range is absent. 



0{L), so that the in- 



a. Small Stokes number differences (0^1) 

In the intermediate asymptotic range r^, <^ r L, 
the particle separation is essentially governed by the first 
term in l|13|) . The quantities of interest are thus described 
by the two-point dynamics associated to a single Stokes 
time - the mean Stokes time. The geometrical inter- 
pretation of this dynamical effect is evident from the 
snapshot of particle positions shown in Fig. |31 Indeed 
at scales larger than the effects of polydispersion are 
negligible. In this range of scales we recover the results 
from the previous Section, namely pi^2ir) ~ Psi^) f^nd 
i'^i,2{r) ~ Ks(r) given by ifTTl) and ifT^ . respectively. 

As displayed in the inset of Fig.|31the situation is quite 
different for r <C <C Now, the time derivative of W 
does not depend on R, but only on the differential accel- 
eration induced by distinct Stokes numbers. Because of 
the resulting relative shift between the two attractors in 
phase space, the velocity difference is independent of the 
separation. As a consequence, impurities with different 
S see each other, at such scales, as a gas of uniformly 
distributed free-streaming particles. This leads to 



Pi,2(r) 



Kl,2(?') 



Csl^l 



fiiS)-d+l 



GsYs\e\ 



r 
L 



L 



(i) 



(15) 
(16) 



These expressions match the intermediate asymptotics 
(|ll|l and H12|) for r = r^,. The constants Cs and Vs are 
related to the two-point motion associated to the mean 
Stokes number S = {Si + S2)/2. 

The presence of a characteristic length scale sepa- 
rating these two regimes is confirmed by numerical ex- 
periments as shown in Fig. ^ 

b. Large Stokes number differences {9 = 0(1)) 

In this case r^, is of the order of the domain size L and 
the intermediate asymptotic regime disappears. There- 
fore, at any scale r <^ L each of the two particles sees the 
other as if uniformly distributed with independent veloc- 
ity. The probability that two particles are at a distance 
r may thus be written 



C 



suS2 n_ 

L \L 



d-l 



(17) 



Regarding the rate of approach, the velocity difference 
can be approximated as the difference between two un- 
correlated velocities. This yields to estimate its typical 
value as (V|^ -|- V|^)^/2^ where Ysj is the characteristic 
velocity of particles with Stokes number Sj. It follows 



Ki,2{r) 



C 



S1.S2 



(v|.+v|.)"= 



(l) 



(18) 
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FIG. 4: Rescaled cumulative probability Pi,2(r) that two par- 
ticles are at a distance smaller than r (lower curve), and cu- 
mulative approaching rate /Ci.2(r) (upper curve) as a function 
of r/r-t. Results are shown for S = 0.18 for different values of 
9 as in the label. Note the collapse of the different curves; the 
dashed lines indicate the scaling behavior below and above r* . 
It is worth also noticing that for r < r* both P and IC scale as 
while for r > r* they scale as r^^^^'^^ and r'''-^-'"'"^, respec- 
tively. As one can see for such a small values of 5* 7 ~ -f 1 . 
Data refer to simulations in the random shear flow llOII with 
T=l, L = l, U = V2S. The curve for Pi, 2 has been shifted 
down for plotting purposes. 



VI. PHENOMENOLOGICAL MODEL FOR THE 
COLLISION KERNEL 

Let us, for the sake of simplicity, focus on polydisperse 
suspensions of particles having the same mass density. 
This assumption implies a one-to-one correspondence be- 
tween the Stokes number and the particle size, namely 
S oa a^. In this framework the simplest phenomenologi- 
cal description of particle collisions in dilute suspensions 
focuses on the total number of collisions per unit time 
between particles of sizes ai and 02 averaged over the 
fluid velocity realizations: 



A/'c(ai, 02) = Q(ai, a2)7Vi A/2, 



(19) 



where A/i is the mean number of particles with size ai in- 
side the domain. The effective collision kernel, Q{ai, 02), 
is given by the average rate at which particles asso- 
ciated to Stokes number Si oc of and S2 oc ar- 
rive at a distance equal to the sum of their radii, i.e. 
Q{ai,a2) = ^1,2(01 + 0'2)- Building quantitative mod- 
els of the collision kernel in turbulent flows is of utmost 
importance for many natural phenomena and industrial 
processes. For instance, we stress that the knowledge 
of Q is crucial to understand the evolution of the droplet 
size distribution in clouds 5^ which is of paramount impor- 
tance for a theoretical comprehension of the rain drops 
formationi^ We believe that the main ingredients needed 
for such a program are, at least at qualitative level, inde- 
pendent on the complexity of the flow. Therefore, here 
we shall concentrate on the simpler case of random flows, 
ignoring some (though important) aspects of turbulent 



flows. Our aim is to develop a semi-quantitative phe- 
nomenological model for (5(01,02) able to capture the 
main ingredients of collisions in suspensions with a broad 
distribution of particle sizes. 

There are two asymptotic regimes where the properties 
of Q are well understood, namely for vanishing inertia 
(S —^ 0) and for infinite response time (S* — *■ 00). As 
stated in the previous Sections, in both limits particles 
distribute uniformly, so that investigating the collision 
kernel reduces to understanding the statistics of velocity 
differences. 

Following Saffman and Turner—, in smooth flows and 
for S* — > the collision kernel can be expressed as 



QsT{ai,a2) = D A ((oi + a2)/£)° 



(20) 



This is obtained by multiplying the geometrical cross sec- 
tion oc (ai -I- 02)''"^ with the typical velocity difference 
between the two particles. The latter is approximated by 
A(ai -I- 02), being A the characteristic fluid velocity gra- 
dient; in turbulent flows the latter is usually estimated 
as A w (e/i/)i/2. 

For 5 ^ 00, it was suggested by Abrahamson^^ that 
the kernel can be obtained by a molecular-chaos type of 
argument, because positions and velocities of particles 
are uncorrelated: 



(ai, 02) 



1? 



in. 



1/2 

S2) (ai 



02) 



(21) 



where Vs ~ Urms/ S^/"^ . 

In the last few years most of theoretical efforts focused 
on the intermediate regime, for which many models and 
predictions have been proposed. For example, in the 
regime S* <C 1 of very small Stokes numbers, collision 
rates are enhanced solely by preferential concentration 
through an increase in the probability of having particles 
at a colliding distance. This was indeed demonstrated by 
means of DNS and theoretical analysis i^ii^iiSiSiiS Thus 
most of the efforts focused in predicting the net effect of 
clusteringi^ This is surely crucial for (close to) monodis- 
perse suspensions of particles with very small sizes, but 
this approach cannot catch collisions between particles 
with different sizes. Moreover, even for same-size par- 
ticles, as soon as the Stokes number reaches small but 
finite values, inertia besides increasing the probability 
to find close particles also enhances particle relative ve- 
locities. This is clear from Fig. ^ where we see that, 
as S approaches the value where clustering is maximal, 
the particle velocity difference is no more proportional 
to the particle separation. This leads to an additional 
increase in the collision rates. Indeed, from comparisons 
with DNS, Wang et al^ argued that most models are 
not able to accurately predict the relative velocity. As to 
collisions between particles having different sizes the situ- 
ation is even more involved. Clearly as soon as the Stokes 
number difference is not negligible the accumulation ef- 
fect induced by inertia need to be understood in terms of 
correlations among particle positions belonging to differ- 
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ent populations (what we called correlations between dif- 
ferent attractors). As far as we know, only few attempts 
in this direction have been considered. Among them we 
mention Zhou et al.-- who investigated by means of DNS, 
bidisperse suspensions in frozen turbulent flows. In their 
study they found a reduction of the accumulation effect 
and an increase in relative velocity. In particular, they 
proposed a simple model for the particle velocities corre- 
lation, in which they prescribed an extremely simplified 
fluid velocity correlation. Though in fairly good agree- 
ment with simulations, their results are rather difficult 
to interpret because spatial correlations of the fluid ve- 
locity are not taken into account. Moreover, their ap- 
proach seems to be justified only for particle response 
time of the order of the correlation time of the large 
scales, Tg. Indeed, they consider in their numerical sim- 
ulations particles with response times up to STg. Also 
Kruis and Kustersli proposed a model for bidisperse col- 
lisions. In their approach they separate, in a somehow ar- 
tificial way (as argued in Ref.[l4ll. two contributions. The 
first named shear mechanism corresponds to the carrier 
fluid velocity shear, that is very similar to the Saffman- 
Turner result. The second, deriving from the acceleration 
induced by the different inertia, was called accelerative 
mechanism. Also this model appear of difficult inter- 
pretation due to the absence of spatial correlation in its 
derivation. Moreover, both mechanisms are considered 
to be effective simultaneously and the dominance of one 
over the other is only due to intensity of the velocity gra- 
dients, while clearly the particle sizes should play a role 
too. 

In the sequel, exploiting the results obtained in the 
previous sections, we propose a phenomcnological model 
of the collision kernel which reduces to the known results 
in the two asymptotics of vanishing and very large Stokes 
number. In the intermediate regime, an interpolation be- 
tween the two asymptotics naturally emerges in terms of 
the exponent 7(5) that, (as seen in Fig.^J out of the two 
limiting regimes, cannot be trivially related to ^. More- 
over, as we shall see, the identification of the crossover 
scale r^, H14|l provides a natural way to discern between 
differential acceleration and shear mechanisms. To be 
used for quantitative predictions an analytical estima- 
tion of 7(5) would be needed. Even for relatively simple 
random flows, this is far from the present capabilities. 
However, this quantity is accessible in simulations and 
may be in principle measured from experiments where 
particle tracking can be done with a high accuracy. 

According to the different cases described in Sec.0 the 
(ai, 02) plane has to be divided into different regions cor- 
responding to the various behaviors of ki,2(?') evaluated 
a,t r — ai+a2- Three regions can be identified (see Fig. [51 
for a sketch): region A where < ai+a-z < L, B that is 
defined by ai -f 02 < r^, < L, and C hy ai + a2 < L < r^,. 
Reminding that oc |aj — a^l/ial -1-02), and the results 
of Sects. IIVI and one has the following behaviors for 
the collision rates. 
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FIG. 5: Sketch of the different regions in the (01,02) plane 
corresponding to different contributions of ki,2 to the effective 
collision rates. 

A- In the gray region defined by 

\ai-a2\ <2L-(4L2-(ai+a2)')'^' , (22) 

the behavior of the mean radial velocity difference is well 
approximated by the two-point motion of particles with 
the mean Stokes number S oc {al+a^)- The collision rate 
thus reduces to that of particles with the same Stokes 
number, namely from Eq. H12|) 

Q(ai,a2) ~ CsVs — j^y(^s)+i — ' 

where the constant Cs depends on the Stokes number 
and on the fluid velocity statistics. 

B- In the white region, the inequality H22|l is fulfilled, yet 
it holds 

|ai-a2| < (2-\/3)(ai+a2). (24) 

The motions of the two particle become uncorrelated and 
from Eq. Hlt)|) one obtains: 

Q(ai,a2)^C5V5 (»2 ^ a2),(5)+i • (25) 

C- In the hatched region where the inequality (|24|l is not 
satisfied, the collision rate is given by 

Q(ai,a2) - Cs^,S2 (Vsi + V5J j2 , (26) 

where Vs^ and are the typical velocities associated 
to particles of size ai and 02, respectively. Of course, 
close to the boundary of the hatched region, the constant 
^Si S2 should be suitably chosen to ensure continuity of 
the collision kernel Q in the plane (01,02). Note that 
this expression is consistent with Abrahamson prediction 
(|21|l . which is expected to hold for large values of the 
Stokes number. 
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FIG. 6; Typical functional shape (bold line) of the collision 
kernel, represented here as a function of 02 through a cut with 
ai = 0.02 fixed. We chose here 7 = 1.8 and the numerical 
factors and constants to fit the order of magnitude obtained in 
our numerical experiments. The effective kernel takes one of 
the different functional forms (I23II . 1251 1 or 12611 (represented 
as dotted lines), depending on whether (ai, 02) is in the region 
A, B or C. 



We now make some comments on the different forms 
taken by the kernel, accordingly to the values of particles 
sizes. The first obvious information is that Eq. (|23|) com- 
prises the two asymptotics of very small (|20|) and very 
large Stokes numbers H21|) when 7 is replaced by its lim- 
iting values. An important observation is that collisions 
between particles with different Stokes numbers are re- 
lated to non-trivial scaling behavior only when their radii 
are rather similar (i.e. in the gray area of Fig. |SJ). In this 
region the rates can be obtained in terms of the dynamics 
for particles with the mean Stokes number. The main in- 
formation is contained in the Stokes-number dependence 
of the exponent 7. 

According to the predictions H23() . (|25|l and H26|l the 
typical shape of the kernel obtained when fixing one of 
the radii and varying the other is represented in Fig. 
The minimal collision rate is obtained for equal-size par- 
ticles; this is consistent with can be understood from the 
symmetry of the kernel under particle exchange. The 
growth of the kernel when 02 3> ai is essentially due to 
the increase of the geometrical cross-section. Note also 
the presence of a maximum when 02 < ai, attained at 
the crossover between (|25|) and (12611 . Note that a similar 
shape has been proposed for the kernel in Ref. 

The numerical results for the collision kernel Q in the 
(01,02) plane are illustrated in Fig. for the Gaussian 
random flow. The one-dimensional cuts represented in 
Fig. |S1 compare favorably with the prediction shown in 
Fig. |B| The random shear flow H10(l displays similar fea- 
tures (not shown). 

In order to quantify the importance of particle inertia 
the ratio Q{a,a)/QsT{a.,a) between the measured ker- 
nel and that obtained from of the Saffman-Turner ap- 
proach H20|l is represented in Fig. |51 To disentangle the 
effects of clustering and densities-velocities correlations. 



FIG. 7: Effective inter-particle collision rate Q{ai,a2) ob- 
tained numerically by considering ghost collisions in the case 
of the time-correlated random flow. The Stokes number S 
and the particle radius a are related by 5 = 2ppa^ /{9pfL^) 
with the choice Pp/pf ~ 4.5 x 10'^. To obtain these rates, 
inter-particle statistics are computed for 100 different values 
of the Stokes number. 
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FIG. 8: Collision rate Q{a\,a2) obtained with the same set- 
tings as for Fig.|7|and represented here for equal-size particles, 
i.e. ai =02, (solid curve, black online) and as a function of 
a2 for two fixed values of ai: a\ — 1.58 x 10"^ (dashed curve, 
blue online), (dotted curve, red online) ai = 3.54 x 10~^. 



the ratio between the kernel obtained when the velocity 
difference is assumed to behave linearly with the separa- 
tion between the two particles, and the Saffman-Turner 
prediction is represented. Note that the two curves coin- 
cide at very small radii: in this regime, the enhancement 
of collision rates is mainly due to clustering effects. Dis- 
crepancies between the two curves appear rather soon, 
before reaching the maximum of clustering. It is clear 
from Fig. |^ that in both cases, there is an increase of 
roughly one order of magnitude in the collision rate of 
inertial particles compared with that of tracers. How- 
ever, the measured values of the kernel differ markedly 
from those obtained when only clustering effects are con- 
sidered. Therefore, away from the two asymptotics it 
is crucial not to consider as independent the effects of 
clustering and enhanced relative velocity. 



11 



^ 11 




0.01 0.02 0.03 0.04 0.05 

Particle radius a 



FIG. 9: Solid curve (black online): collision rate Q{a,a) for 
equal-size particles, normalized by that obtained when ne- 
glecting particle inertia, i.e. the Saffman and Turner result 
Eq. (1201 . Dashed curve (blue online): collision rate obtained 
when neglecting correlations between the velocity difference 
and the density (i.e. assuming that the first is just propor- 
tional to the separation between the particles), normalized in 
the same way. These curves were obtained numerically with 
the same setting as for Fig.Q the Saffman-Turner kernel was 
evaluated by following tracer particles advected by the same 
flow. 

VII. CONCLUDING REMARKS 

We have identified the main mechanisms leading to 
the enhancement of clustering and coUisions induced by 
inertia in dilute suspensions of heavy particles. A La- 
grangian method based on ghost collisions has been used 
for numerical investigations. In agreement with previ- 
ous studies we found that clustering is maximal when 
the Stokes number is order unity. However, this is not 
the only mechanism enhancing collisions: large veloc- 
ity differences can occur at very small separations and 
for finite values of the Stokes number. This results in 
non-trivial scaling properties of the rate at which parti- 
cles come close. A phenomenological model for the colli- 
sion kernel based on these ingredients has been proposed. 
Our results highlight the importance of accounting for 
the full position- velocity phase-space dynamics, particu- 
larly for polydispcrse solutions. It is important to stress 
that, away from the small Stokes asymptotics, multipha- 
sic approaches )^L2^ based on a continuum description of 
suspensions in position space, may fail to catch these ef- 
fects. 

We conclude by discussing some aspects concerning im- 



purities in turbulent flows. 

We expect our results to be relevant to clustering and 
collisions at dissipative scales. Turbulence will certainly 
affect both the values of the constants and the scaling 
exponents, and induce a non-trivial Reynolds number de- 
pendence of the collision rates. Within the framework of 
model flows it would be interesting to extend the present 
study to random flows with non-Gaussian statistics to ac- 
count for intermittency of the velocity gradients. More- 
over, in actual turbulent flows the Kolmogorov time-scale 
is a dynamical variable which has important fluctuations. 
Therefore, the Stokes number is also a random variable 
along the particle trajectories. This provides a further 
complication in developing models accounting for clus- 
tering and collisions of particle suspensions, which needs 
to be investigated. 

Another important issue concerns clustering at 
inertial-range scales, where the velocity field is not dif- 
ferentiable. Experiments indeed show that preferential 
concentration appears also at those scalesi Non-trivial 
clustering properties have been observed numerically also 
in the inverse-cascade range of two-dimensional turbulent 
flows, namely the formation of holes in the distribution of 
particlesi^ It is important to remark that clustering at 
the inertial scales may influence the probability for two 
particles to arrive below the dissipative scale and thus the 
collision rates. In the inertial range the dynamics of the 
fluid is close to Kolmogorov 1941 theory, i.e. the velocity 
fleld is Holder continuous with exponent 1/3. As a con- 
sequence, tracers separate explosively giving rise to the 
celebrated Richardson's i'^/^ law. For inertial particles, 
one needs to understand the competition between explo- 
sive separation and clustering due to dissipative dynam- 
ics. In this direction, it may be useful to further extend 
to inertial particles recent models and techniques devel- 
oped in the framework of passive scalars (for two recent 
reviews see Refs. liollill) . 
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